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We report on our study of QCD thermodynamics with 2+1 flavors of dynamical quarks. In this 
proceeding we present several thermodynamic quantities and our recent calculation of the critical 
temperature. In order to investigate the thermodynamic properties of QCD near the continuum 
limit we adopt improved staggered (p4) quarks coupled with tree-level Symanzik improved glue 
on N t = 4 and 6 lattices. The simulations are performed with a physical value of the strange quark 
mass and light quark masses which are in the range of m q /m s = 0.05 — 0.4. The lightest quark 
mass corresponds to a pion mass of about 150 MeV. 
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1. Introduction and Lattice Setup 

The calculation of QCD thermodynamics from first principle is important for various research 
areas such as Heavy Ion Phenomenology, Cosmology and Astrophysics. Lattice QCD enables us 
to carry out such calculations. Especially for HIC phenomenology it is mandatory to improve 
estimates on some basic thermodynamic quantities which have been obtained in previous lattice 
calculations. One of the main subjects in our project is the accurate determination of the criti- 
cal temperature T c , whose uncertainty, for example, strongly affects the critical energy density e c 
because of its T c dependence, e c ~ T c 4 . 

Since thermodynamics of lattice QCD requires huge computational resources, it is difficult 
to perform an ideal simulation. Recent studies tell us that quark masses and the number of flavors 
strongly affect thermodynamic quantities Reliable continuum extrapolations are of tremendous 
importance as well [0]. Therefore, it is our aim to study QCD thermodynamics with almost real- 
istic quark masses on the QCDOC machine at Brookhaven National Laboratory and the APEnext 
machine at Bielefeld University. The calculation is performed with Nf = 2 + 1, which means 2 
degenerate light quarks and one heavier quark on lattices with N t = 4 and 6. The lightest quark 
masses of our simulation yields a pion mass of about 150 MeV and a kaon mass of about 500 MeV. 

For such calculations we adopt the p4fat3 quark action, which is an improved Staggered quark 
action j3|], with a tree-level improved Symanzik gauge action. By using the p4fat3 action, the free 
quark dispersion relation has the continuum form up to 0(p 4 ), and the taste symmetry breaking is 
suppressed by a 3-link fattening term. The action also improves bulk thermodynamical quantities 
in the high temperature limit [Qj. The improvements are essential to control the continuum extrap- 
olation on rather coarse lattices, i.e. N t = 4 and 6. The gauge ensembles are generated by an exact 
RHMC algorithm [§]. 

As a status report of the project, in this proceeding, we present several thermodynamic quanti- 
ties, which are order parameters and their susceptibilities, the static quark potential, and the spatial 
string tension. In the last section we discuss the critical temperature at the physical point. The 
details of the critical temperature calculation are given in our recent paper [g]. 

2. Order Parameters and Susceptibilities 

To investigate the QCD critical temperature and phase diagram, order parameters of the QCD 
transition are indispensable. In the chiral limit the chiral condensate (yy) is the order parameter 
for the spontaneous chiral symmetry breaking of QCD. On the other hand in the heavy quark limit 
the Polyakov loop (L) is the order parameter of the deconfinement phase transition. For finite quark 
masses, these observables remain good indicators for the (pseudo) critical point. Especially their 
susceptibilities are useful to determine the critical coupling /3 C in numerical simulations. 

Figure [I] shows the susceptibilities of the light quark chiral condensate 1 on 8 3 x 4 and 16 3 x 4 
lattices with various quark masses. The peak positions of the susceptibilities define the point of the 
most drastic change of each order parameter, i.e. the (pseudo) critical point of the QCD transition. 
The results are interpolated in the coupling /3 by using the multi-histogram re-weighting technique 
§• 

'We always use here the disconnected part of the chiral susceptibility. 
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Figure 1: The susceptibility of the light quark chiral condensate on 8 3 x 4 and 16 3 x 4 lattices. Each panels 
correspond to results with thi / m s — 0.05 , 0. 1 , 0.2 and 0.4 respectively. The lines are calculated by the multi- 
histogram re-weighting technique. 
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Figure 2: The difference of gauge couplings at the location of peaks in the Polyakov loop and the chiral 
susceptibilities, f} L — 



The strength of the transition decreases with increasing quark masses, this is reflected in the 
decreasing peak height of the chiral susceptibilities. We calculate these susceptibilities on lattices 
with aspect ratios of N s /N t = 2 and 4. Since we see a rather small volume dependence the results 
suggest that the transition is in fact not a true phase transition in the thermodynamic sense but a 
rapid crossover. 

Figure ^| (left) shows the difference in the peak position of the chiral and Polyakov loop sus- 
ceptibilities. The differences are small and almost identical at the aspect ratio of N s /N t = 4, which 
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Figure 3: Heavy quark potentials scaled by ro for various cut-offs and quark masses (left) and a global fit of 
ro/a with the RG inspired ansatz Eq.(3.3) (right). 



indicates the chiral and deconfinement phase transition occur at almost the same temperature. The 
discrepancy between the peak positions shrinks with increasing volume. 

3. Scale Setting and the Heavy Quark Potential 

The lattice scale is determined by the heavy quark potential V(r) which is extracted from 
Wilson loops. The Wilson loop expectation values are calculated on 16 3 x 32 lattices with APE 
smearing in spatial direction. The spatial path in a loop is determined by the Bresenham algorithm 
[S]. We calculate the string tension, a and Sommer scale ro, which is defined [^] as the distance 
where the corresponding force of the static quark potential matches a certain value suggested by 
phenomenology: r 2 ^ \ r = ro = 1.65. To remove short range lattice artifacts we use the improved 
distance, r !mp , which is defined as 

1 f d 3 k e ikr 



(3.1) 



Anr imp J (27r) 3 4I i -(sin 2 |) + isin 4 |' 

In our lattice setup, we find almost no mass and cutoff dependence in the potential scaled by ro 

at ,/V T = 4 and 6 (Fig.2(left)). As discussed in previous studies we also find no string breaking 

effects even at large r. To estimate systematic uncertainties of the potential fit, we performed 

several types of fits, e.g. different fit-ranges in r and fit-forms (3 & 4 params. fits), 

oc oc oc oc 
V{r)=C+ + Cr imp , V{r)=C+- + Or + d{ ). (3.2) 

l~imp r I~j m p V 

The differences in the mean values of the fits are taken into account as a systematic uncertainty of 
the scale setting. 

We have determined the scale parameter r in units of the lattice spacing for 9 different pa- 
rameter sets. This allows to interpolate between different values of the gauge coupling and quark 
masses. We use a renormalization group inspired ansatz [[H]] which takes into account the quark 
mass dependence of ro/a [[y]] and which approaches, in the weak coupling limit, the 2-loop /3- 
function for three massless flavors, 

(r /a)- 1 =R(fi)(l +Ba 2 (p) + Co 1 (p)) e A ( 2A i+^)+ D . (3.3) 
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Figure 4: Temperature dependence of the spatial string tension for Nf = and Nf = 2 + 1. Both dotted lines 
are fits with Eq. (4.1). Note, that the scale on the horizontal axis is logarithmic. 



Here /?(/$) denotes the 2-loop j6-function and = /?(/$)//?(/$) with j8 = 3.4 chosen as an arbi- 
trary normalization point. 



4. Spatial String Tension 



Let us now discuss the calculation of the spatial string tension which is important to verify 
the theoretical concept of dimensional reduction at high temperatures. The spatial string tension 
is extracted from the spatial static quark "potential" (from spatial Wilson loops). We use the same 
analysis technique as for the usual (temporal) static quark potential. 

At high temperature, the spatial string tension o s (T) is expected to behave like 



x/a^T) = cg 2 (T)T. 

Here g 2 (T) is the temperature dependent coupling constant from the 2-loop RG equation, 



(4.1) 



-\T) =2Mn-^ + ^ln(21n-^ 
A ff b Q V A a 



(4.2) 



If dimensional reduction works, the parameter"c" should be equal to the 3-dimensional string ten- 
sion and should be flavor independent. 

Our 2+1 flavor result yields c = 0.587(41) and A a /T c = 0.114(27), obtained by a fit with 
Eq.(|Q|). On the other hand, we plot in Fig. || also the quenched result [ |l2| ] which gives c = 
0.566(13) and A a /T c = 0.104(9). We thus find that the parameter "c" is - within statistical errors 
- independent on the number of dynamical flavors and that dimensional reduction works well even 
for T = 2T C . This analysis can and will be refined in the future by taking into account higher order 



corrections to Eq.(4.1) [13]. 
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Figure 5: T c ro (left) and T c /^/a (right) as a function of m ps ro on lattices with temporal extent N T = 4, 
m s = 0.065 (squares) and m s = 0.1 (triangles) as well as for N z = 6, m s = 0.04 (circles). Thin error bars 
represent the statistical and systematic error on ro/a and \faa. The broad error bar combines this error with 
the error on j3 c . The vertical line shows the location of the physical value m ps rQ = 0.321(5) and its width 
represents the error on r$. The three parallel lines show results of fits based on Eq. (5.1) with d = 1.08 for 
N T = 4, 6 and Nz — ► °° (top to bottom). 



5. The transition temperature 

Finally I discuss the transition temperature in QCD, which is one of the most important sub- 
jects in our project. In sect.|| we have determined the critical /3 at each N T and for several quark 
masses. At these couplings we performed zero temperature calculations of the static quark poten- 
tial. The scale settings at each critical j3 provide each values for the critical temperature. We thus 
can extrapolate T c to the chiral limit as well as to the physical point by using a scaling ansatz. In 
the zero temperature calculations the actual j8 is sometimes slightly different from our final result 
on the critical j8. The differences are corrected by using Eq j3.3| . Furthermore a systematic error in 
the critical /3, e.g. jS/, — /3/ at N T = 6 in Fig.||, is also taken into account in the critical temperature. 

In Fig. |5] we plot the critical temperature in unit of the Sommer scale (left) and the string 
tension (right) as function of the pion mass mps (also in units of the Sommer scale). Thin error 
bars represent the statistical and systematic error on ro/a and \[ca. The broad error bar combines 
this error with the error on j8 f . We perform a combined chiral and continuum extrapolation of T c 
by using the ansatz, e.g. in unit of the Sommer scale, 

T c r = [T c r ]^f +A{m n r ) d + B/N^, (5. 1) 

where A and B are free fit parameters. If the QCD transition is second order in the chiral limit 
the transition temperature is expected to depend on the quark mass as m\^ S , which corresponds 



to d ~ 1 .08 in Eq. |5.1| , by the fact that one expects a critical point in the chiral limit which is in 
the (?(4)-universality class. If, on the other hand, the transition becomes first order for small quark 
masses, the transition temperature will depend linearly on the quark mass, i.e. d = 2. Using the fit 
form ansatz |5j] we can determine the transition temperature at the physical point, which is defined 
as m ps ro = 0.321(5), 

T c r = 0.457 (7)^ 2 , T c /^ = 0.408(8)1?, (5.2) 
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where the central value is given for fits with d = 1.08 and the lower and upper systematic error 
correspond to d = 1 and d = 2, respectively. 

In order to convert T c to physical units, we use the scale parameter, ro = 0.469(7) fm, deduced 
from the bottomonium level splitting using NRQCD Jl4|]. Finally we obtain the transition tem- 
perature in QCD at the physical point, T c = 192(7) (4) MeV, where the statistical error includes 
the errors given in Eq. ||j2| as well as the uncertainty in the value of ro and the second error re- 
flects our estimate of a remaining systematic error on the extrapolation to the continuum limit. 
The value of the critical temperature obtained here is about 10% larger than the frequently quoted 
value ~ 175 MeV. We note that this larger value mainly results from the value for ro used in our 
conversion to physical scales. 

The analysis presented here leads to a value for the critical temperature with about 5% statisti- 
cal and systematic errors. It clearly is desirable to confirm our estimate of the remaining systematic 
errors through an additional calculation on an even finer lattice. Furthermore, it is desirable to 
verify this result through calculations using other T = scales and to explore other discretization 
schemes for the fermion sector of QCD and to also obtain a reliable independent scale setting for 
the transition temperature from an observable not related to properties of the static quark potential. 
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